Sparse deconvolution (an MM algorithm)

Deconvolution refers to the problem of estimating the unknown input to an LTI system when the output signal and system response are known. In practice, the available output signal is also noisy. For some systems, the deconvolution problem is quite straight forward; however, for systems that are non-invertible or nearly non-invertible (e.g. narrow-band or with frequency response nulls), the problem is more difficult. The use of an exact inverse system can greatly amplify the noise rendering the result useless. In such cases, it is important to utilize prior knowledge regarding the input signal so as to obtain a more accurate estimate of the input signal, even when the system is nearly non-invertible and the observed output signal is noisy. In some applications of deconvolution, it is known that the input signal is sparse (i.e. a spike train, etc.) or approximately sparse. Applications of ‘sparse deconvolution’ include geophysics, ultrasonic non-destructive evaluation, speech processing, and astronomy [11]. One approach to sparse deconvolution involves the minimization of a cost function defined in terms of the `1 norm [4, 6, 16]. The minimization of cost functions defined in terms of the `1 norm is useful not just for deconvolution, but for sparse signal processing more generally. Indeed, since its early application in geophysics, the `1 norm and sparsity have become important tools in signal processing [3]. The tutorial [14] compares least squares and `1 norm solutions for several signal processing problems, illustrating the advantages of a sparse signal model (when valid). This tutorial aims to illustrate some of the principles and algorithms of sparse signal processing, by way of considering the sparse deconvolution problem. A computationally efficient iterative algorithm for sparse deconvolution is derived using the majorization-minimization (MM) optimization method. The MM method is a simple, yet effective and widely applicable, method that replaces a difficult minimization problem with a sequence of simpler ones [8]. Other algorithms, developed for general `1 norm minimization, can also be used here [5,13,17]. However, the MM-derived algorithm takes advantage of the banded structure of the matrices arising in the sparse deconvolution problem. The resulting algorithm uses fast solvers for banded linear systems [1], [12, Sect 2.4]. The conditions that characterize the optimal solution are described and illustrated in Sec. 3. With these simple conditions, the optimality of the result computed by a numerical algorithm can